Population pharmacokinetics and limited sampling strategy for therapeutic drug monitoring of mycophenolate mofetil in Japanese patients with lupus nephritis

Background Mycophenolate mofetil (MMF), a prodrug of the immunosuppressive agent mycophenolic acid (MPA), is difficult to administer because of the pharmacokinetic complexity of MPA. Although dosage adjustment according to the 12-h area under the concentration–time curve (AUC0-12) is thought to be desirable, multiple blood samplings for AUC calculation may pose a clinical challenge. A limited sampling strategy (LSS) would provide a solution; however, little is known about MPA pharmacokinetics in lupus nephritis patients, especially in those with Asian backgrounds, or few, if any, LSSs are reported for them. Methods Thirty-four adult Japanese patients receiving MMF for lupus nephritis were examined retrospectively. MPA pharmacokinetics were investigated, and a PPK model was developed using Phoenix® NLME™ software. Single and double blood sampling strategies from Bayesian estimation using the PPK model and from multiple linear regression were compared. Tolerability was also evaluated. Results In the pharmacokinetic analysis, renal function and serum albumin had significant effects on dose-normalized AUC0-12; and serum albumin, concomitant proton pump inhibitor (PPI) and iron/magnesium oxide did on dose-normalized maximum concentration. As a PPK model, a two-compartment model was developed with a transit absorption model and first-order elimination, in which creatinine clearance and serum albumin were covariates for MPA clearance. The double sampling strategy at 1 and 4 h by multiple linear regression showed the best agreement with the observed AUC0-12 (r2 = 0.885). Of the single sampling strategies, the one at 6 h by Bayesian estimation performed best (r2 = 0.769). The tolerability evaluation showed that correlations were suggested for gastrointestinal involvement. Conclusions The present study developed the first PPK model of MPA for Japanese lupus nephritis patients. As for LSSs, a double sampling strategy at 1 and 4 h by multiple linear regression would work best; when only a single blood sampling is allowed, a strategy at 6 h by Bayesian estimation using the PPK model developed in this study would be best. The LSSs good enough for clinical use may facilitate safer, more effective, and individualized therapy. Supplementary Information The online version contains supplementary material available at 10.1186/s40780-022-00271-w.


Background
Mycophenolate mofetil (MMF) is a prodrug of the immunosuppressive agent mycophenolic acid (MPA). MMF is widely used in kidney and other organ transplantations and is also recommended as induction and maintenance therapy for lupus nephritis (LN), a renal manifestation of systemic lupus erythematosus [1]. Orally administered MMF is almost completely absorbed by the gastrointestinal tract and rapidly hydrolyzed to its active form MPA, which inhibits inosine 5'-monophosphate dehydrogenase (IMPDH), an enzyme involved in the de novo synthesis of guanosine in lymphocytes [2]. MPA is then primarily metabolized into the pharmacologically inactive mycophenolic acid glucuronide (MPAG) in the liver [3] and thereafter hydrolyzed back to MPA during enterohepatic recirculation (EHC), resulting in a second peak in MPA concentration [4].
In administering MMF, the pharmacokinetics of MPA should be taken into consideration. MPA has a narrow therapeutic window and exhibits wide inter-and intraindividual pharmacokinetic variability. The pharmacokinetics are affected by numerous factors, including renal function, liver function, serum albumin levels, concomitant drugs, and ethnic background [4][5][6][7][8][9][10][11]. Dosage reduction or drug withdrawal is sometimes necessary when adverse effects arise, such as infections or gastrointestinal symptoms [12][13][14][15]. Although a few studies have reported correlations between pharmacokinetic parameters and adverse effects [16][17][18], much remains unknown. Because of the administrative difficulties, adjusted dosage is desired rather than fixed dosage.
The 12-h area under the concentration-time curve (AUC 0-12 ) for MPA reportedly correlates well with clinical outcomes [19,20] and is considered to represent a valuable tool for adjusting MMF dosages. Although trough concentration could provide an alternative to AUC in general, such single-point measurements would not be appropriate in the case of MMF due to the second peak in MPA concentration [4,21]. In patients receiving renal transplantation, a target range of 30-60 μg•h/ mL in AUC 0-12 has been proposed in general, and therapeutic drug monitoring (TDM) using AUC 0-12 has been suggested as a strategy for personalized treatment [22]. One clinical challenge with AUC 0-12 is that calculation requires multiple blood samples. In MMF, 8-10 blood samples within 12 h would be suggested for ideal calculation, but in reality, blood collection at such a high frequency is not only costly and time-consuming, but also places a huge burden on patients. As a solution, a limited sampling strategy (LSS), which estimates AUC 0-12 from a small number of samples, may be developed. For transplant patients receiving MMF, LSSs, varying from patient group to patient group, are commonly applied in TDM using AUC 0-12 [23].
For Asian LN patients treated with MMF, however, the pharmacokinetic properties of MPA remain largely unexplored and few, if any, population pharmacokinetic (PPK) models, which may help in the development of LSS, appear to have been built. PPK models have been reported for transplant patients with various backgrounds [24][25][26], and a few have been described for autoimmune patients of different ethnicities [10,27]; but these may not be applicable to Asian LN patients because of potential differences in MPA pharmacokinetics. While previous studies have reported that a range of 30-45 μg•h/mL in AUC 0-12 is associated with good clinical outcomes [28], an LSS appears less likely to be available for Asian LN patients.
This study aimed to develop an LSS for Japanese LN patients receiving MMF, in the hope that the LSS would help facilitate safer, more effective, and individualized therapy for such patients. For this purpose, MPA pharmacokinetics after oral administration of MMF were investigated first, and then a PPK model was constructed to describe the pharmacokinetics. From the resulting model, an LSS was developed to estimate AUC 0-12 from single or double blood sampling. For comparison, two distinct approaches were employed in developing an LSS: maximum a posteriori Bayesian estimation (MAP-BE) and multiple linear regression (MLR).

Patients and data collection
Participants comprised all Japanese patients ≥ 18 years old who received MMF (CellCept ® ; Chugai Pharmaceutical Group Co., Tokyo, Japan) for treatment of LN at Aichi Medical University Hospital (Nagakute, Japan) between March 2015 and June 2022 and for whom data were available on multiple blood samples necessary to make a pharmacokinetic profile.

Pharmacokinetic analysis
The parameters used for pharmacokinetic analysis were maximum MPA concentration within a dosing interval of 12 h (Cmax), time to reach maximum MPA concentration (Tmax), and AUC 0-12 . Cmax and Tmax were derived from the MPA concentration-versus-time profiles of patients. AUC 0-12 was calculated from MPA concentrations at the 9 time points of C 0 to C 12 , using the linear trapezoidal method. Dose-normalized Cmax and dosenormalized AUC 0-12 were also calculated. Correlations were investigated between AUC 0-12 and MMF dose, as well as between AUC 0-12 and C 0 MPA concentration, using the Spearman rank correlation test, where AUC 0-12 was regarded as a dependent variable and the dose or concentration as an independent variable.
From a clinical point of view, we performed group comparisons to investigate the effect of biochemical factors on the pharmacokinetic parameters of MPA, as the previous study [31] did. The patients were divided into two groups by eGFR (< 81 mL/min versus ≥ 81 mL/ min), and also by serum albumin level (< 3.5 g/dL versus ≥ 3.5 g/dL). Similarly, for the effect of concomitant drugs, patients were split into two groups by the presence or absence of each of the following drug administrations: tacrolimus, PPI, and iron/magnesium oxide. In all group comparisons, differences in Tmax, dose-normalized Cmax, and dose-normalized AUC 0-12 were examined using the Mann-Whitney U test.
All statistical analyses of pharmacokinetics were performed using easy R (EZR) version 1.54 [32].

PPK model development and validation
All analyses for PPK model development and validation, as well as subsequent LSS development and validation, were performed using Phoenix ® NLME ™ software (version 8.1; Pharsight, Mountain View, CA, USA) and EZR version 1.54 [32]. For the estimation algorithm, the firstorder conditional estimation-extended least-squares estimation method was employed. For the overall procedure, Bonate [33], Owen and Fiedler-Kelly [34], Gabrielsson and Weiner [35], and Guidance for Industry: Population Pharmacokinetics [36] were referred to.
The PPK model was developed in a step-by-step manner. First, in the search of a base model, one-, two-and three-compartment models with first-order elimination were tested. These models were parameterized in terms of absorption rate constant (Ka), central volume of distribution (V 1 ), central compartment clearance (CL), peripheral volume of distribution (V 2 ), and intercompartmental clearance (Q). Because absolute bioavailability (F) was not to be assessed in cases of oral administration, V 1 /F, CL/F, V 2 /F, and Q/F were considered to correspond to V 1 , CL, V 2 , Q, respectively. Using the resulting base model, the intra-individual variability of MPA concentration was examined by comparing additive, proportional, and mixed (additive plus proportional) residual error models. Inter-individual variability was examined using the exponential error model.
As the next step, covariate models were developed and the following variables were tested as covariates: age, sex, body weight, eGFR, CLcr, serum albumin, UPC, ALT, AST, TBIL, CRP, dose of concomitant prednisolone, and concomitant drugs of PPI and iron/magnesium oxide. In addition, an EHC model was also examined in reference to a previous study [27]. Finally, whether the MPA absorption process was well described was examined by comparing models with and without lag time. For the same purpose, transit models with different numbers of compartments (1 to 7) were also compared.
In these steps of model development, model selection decisions were made based primarily on the likelihood ratio test. In the test, objective function values (OFVs) in nested models were compared using the χ 2 test. If the difference between models resulted in OFV greater than the critical value (i.e., 6.63 units when models differ by 1 degree of freedom), the more complex model was considered to be significantly better (p < 0.01). Subsequently, relative standard errors (RSE) were calculated for the evaluation of adequacy, and the condition number (i.e. the ratio of the absolute highest and lowest eigenvalues) for overparameterization. In both calculations, a smaller value was considered to be better. As another way to see adequacy, goodness-of-fit plots were performed. Goodness-of-fit plots were examined for observed concentrations versus individual predicted concentrations, observed concentrations versus population-predicted concentrations, conditional weighted residuals versus time since the dose, and conditional weighted residuals versus population-predicted concentrations [37,38]. For validation, models in the development process were examined using the prediction-corrected visual predictive check (pcVPC, 1000 simulations) [39]. The precision of parameter estimates was assessed using the bootstrap method (1000 bootstrap samples).

LSS development and validation
From a clinically practical point of view, we developed single and double blood sampling strategies. In the development, we employed two distinct methods: an MAP-BE approach using the PPK model and an MLR approach [40]. In the MAP-BE approach, from one or two observed MPA concentrations, the remaining values from C 0 to C 12 were estimated based on the newly developed PPK model. From the observed and estimated MPA concentrations, the predicted value of AUC 0-12 (AUC predicted ) was calculated. In the MLR approach, multiple linear regression analysis was performed between one or two observed values of MPA concentration (independent variables) and AUC observed (dependent variable), which was calculated in the same way as AUC 0-12 . Using the resulting MLR equation, AUC predicted was calculated from the one or two observed values.
For validation, adjusted coefficients of determination (r 2 ) were used to examine the regression level of AUC predicted and AUC observed . Accuracy and precision of prediction were evaluated at each time point from C 0 to C 12 , using mean prediction error (MPE), 95% confidence interval (CI) of MPE, and root mean squared percentage error (RMSPE). Acceptance ranges were set as ≤ ± 10% in MPE and ≤ 25% in RMSPE [40]; CI of MPE was considered to be acceptable when it crossed the value of 0. In the equations below, N represents the number of pairs of AUC predicted and AUC observed .

Tolerability evaluation
As major adverse effects, gastrointestinal symptoms and infections were examined. Using the Mann-Whitney U test, patients with and without gastrointestinal involvement were compared in terms of dose, dose/body weight, C 0 , Cmax, AUC 0-12 , AUC 0-12 /dose, and AUC 0-0.5 ; and so were patients with and without infection. If the correlation between a parameter and either of the major adverse occurrences was suggested, receiver operating characteristic (ROC) curve analysis was performed for the parameter, to determine the optimal cutoff. Statistical analyses were performed using EZR version 1.54 [32].

Patients
Thirty-four Japanese ≥ 18 years old who received MMF for treating LN were enrolled. They received MMF at doses of 250-1,000 mg every 12 h twice a day, as prescribed in the package insert. Demographic and clinical characteristics of the patients are shown in Table 1.
In 33 patients, data for blood samples taken at all 9 time points (C 0 to C 12 ) were available, making full pharmacokinetic profiles. In 1 patient, data for only 5 time points (C 0 , C 1 , C 2 , C 3 , and C 4 ) were obtained. The data set contained 302 MPA plasma concentrations obtained from 34 patients, which were used for PPK model construction. For pharmacokinetic analysis and LSS development, full profiles of the 33 patients were used.
The effects of concomitant administration of PPI nor iron/magnesium oxide were not observed in dose-normalized AUC 0-12 (Fig. 2C D), but in dose-normalized Cmax. It was significantly lower in patients with concomitant PPI (n = 21) than in those without it (n = 12) (Fig. 2G); and likewise in patients with iron/magnesium oxide (n = 5), compared to those without it (n = 28) (Fig. 2H). On Tmax, concomitant iron/magnesium oxide showed an effect: patients receiving iron/magnesium oxide (n = 5) displayed a significantly increased Tmax compared to those without such treatment (n = 28).

PPK model development and validation
The structure of the final PPK model that adequately described MPA pharmacokinetics of the 34 patients is given in Fig. 3. As a base model, a two-compartment linear model with first-order elimination was employed.
Intra-individual variabilities were best described by a proportional residual error model. Then models with different covariant patterns were developed. Of these models, the one including all the variables significant in the pharmacokinetic analysis, i.e. CLcr and serum albumin   (Table 2 and Additional file 2). The model with CLcr and serum albumin alone as covariates was thus adopted for the final model. In the next step, transit compartments were added one at a time, up to seven, and each model was examined. The model with lag time and the EHC model were also examined. The model with 6 compartments significantly improved the model fit and better described the absorption process. The model with lag time or the EHC model did not improve the model fit.
All goodness-of-fit plots showed the high predictive performance of the final model, with no systematic deviations observed (Fig. 1B-E). The relative errors were comparable across all concentration ranges (Fig. 1E,  Additional file 3). The prediction-corrected visual predictive check is shown in Fig. 1F. The observed median and 5th and 95th percentiles were located inside the 90% CI of the simulated data. Parameter estimates by the bootstrap method closely matched the means of corresponding parameter estimates from the final model, confirming the stability of the final model ( Table 2). Table 3 shows the best four predictive performances for each of single and double blood sampling strategies under the MAP-BE approach (using Bayesian estimation), and the MLR approach (using multiple linear regression analysis). For both approaches, double sampling strategies were more accurate than single ones. In the MAP-BE approach, the MPE and RMSPE values were all within the clinically acceptable ranges of ≤ 10% and ≤ 25%, respectively; and the CIs of MPE crossed the value of 0, except for single sampling at C 0 , single sampling at C 4 , double sampling at C 1 and C 3 , and double sampling at C 1 and C 4 . In the MLR approach, the RMSPE values of the single sampling strategies were all larger than the acceptable range, while those of the double sampling strategies were all within the acceptable range. The MPE values were all within the acceptable range and the CIs of MPE all crossed the value of 0. The best agreement with AUC observed was achieved by the double blood sampling strategy of C 1 and C 4 by the MLR approach, where r 2 was 0.885, MPE was 2.50% (95%CI, -2.74-7.73), and RMSPE was 14.75%. Of the single sampling strategies, C 6 sampling strategy by the MAP-BE approach displayed the best agreement with the AUC 0-12 , where r 2 was 0.769, MPE was -0.72% (95%CI, -7.51-6.08), and RMSPE was 18.89%. The estimation accuracy was found to be good at near the target value of AUC 0-12 (30-45 μg•h/mL) (Fig. 4).

Tolerability evaluation
The results for major adverse occurrences (gastrointestinal involvement and infection) are shown in Table 4. Gastrointestinal involvement occurred in 10 of the 33 patients, and diarrhea, vomiting, and epigastric pain were observed. Infections occurred in 14 patients, and upper respiratory infection, fever, pharyngitis, bronchitis, oral herpes zoster, Malassezia, and Aspergillus were observed. No case of cytomegalovirus infection was observed. For gastrointestinal involvement, correlations were suggested for three parameters: AUC 0-0.5 , Cmax and AUC 0-12 were significantly higher in patients with gastrointestinal involvement (n = 10) than those without gastrointestinal involvement (n = 23) (p < 0.01, p < 0.05 and p < 0.05, respectively). No correlation was suggested between infection and any of the examined parameters.
ROC curve analysis was performed for each of the three parameters for which correlation with gastrointestinal involvement was suggested. The analysis indicated the optimal cutoff for each parameter: 3.880 μg•h/mL for AUC 0-0.5 , 24.360 μg/mL for Cmax, and 47.597 μg•h/mL for AUC 0-12 (Additional file 4A-C).

Discussion
This study retrospectively investigated the electric medical records of 34 Japanese LN patients receiving MMF. MPA pharmacokinetic analysis revealed that renal function, serum albumin, concomitant PPI, and iron/ magnesium oxide had significant effects (Fig. 2). As a PPK model, a two-compartment model with first-order elimination was developed in which the absorption process was best described by fixed 6 transit compartments (Fig. 3), and intra-individual variabilities were best described by a proportional residual error model. In the model, CLcr and serum albumin were significant covariates for CL (Table 2). For an LSS to be developed, single and double blood sampling strategies were examined by MAP-BE and MLR approaches. The results showed that the double sampling strategy of C 1 and C 4 by MLR best agreed with AUC observed (r 2 = 0.885). Of the single sampling strategies, the one at C 6 by MAP-BE agreed best (r 2 = 0.769) ( Table 3). MPA pharmacokinetics are reportedly affected by renal function. Reduced renal function can affect the MPA pharmacokinetics, leading to reduced elimination of and thereby higher exposure to MPA [4]. Because around 97% of MPA binds to albumin, hypoalbuminemia increases the free fraction of MPA, resulting in faster clearance of and thereby reduced exposure to MPA. In fact, low serum albumin concentrations have been reported to increase CL [4,6]. In this study too, covariate analysis identified the effects of CLcr and serum albumin on CL. Although MPA pharmacokinetics could differ among populations [40], the estimated CL in our study was 13.15 L/h, comparable with those in the previous studies of Thai (14.5 L/h) [41] and Mexican (15.4 L/h) [10] patients with LN. Taken together, the MPA pharmacokinetics of the population in this study were similar in many ways to those in the previous studies: dose-normalized AUC 0-12 increased significantly when renal function was damaged, but decreased significantly when serum albumin concentrations were low [40,42]; and dose-normalized Cmax decreased significantly with concomitant PPI and iron/magnesium oxide (Fig. 2) [7,43,44]. These effects were taken into consideration in the first PPK model, in which creatinine clearance and serum albumin were covariates for MPA clearance and concomitant PPI and iron/magnesium oxide were covariates for V 1 (Additional files 2, 5). Absorption of MPA increases at low intragastric pH as the dissolubility of MMF increases [7,45]. Although PPI was found to be a possible covariate of V 1 , histamine H2-receptor antagonist, which restrains gastric acid secretion, was not. A plausible speculation is that concurrent use of histamine H2-receptor antagonists does not affect pharmacokinetics overly much because these agents wear off as time passes, whereas PPIs impede irreversibly and provide lasting effects [46]. Besides PPI, iron/magnesium oxide decrease absorption of MPA as it forms a chelate with MMF [5,47]. These facts may explain why concomitant PPI and iron/magnesium oxide were identified as possible covariates on V 1 .
PPI and iron/magnesium oxide, which are not always used for treating LN, were coadministered in 22 (64%) and 5 (14.7%) out of 34 patients, respectively, in this study. They were found to have significant effects in the pharmacokinetic analysis, and the PPK model including them, in addition to CLcr and serum albumin, as covariates significantly decreased OFV. Yet, suspecting overparameterization and considering evaluation adequacy and precision, we decided to exclude them from covariates for the final model. For confirmation of the decision, we examined LSSs from Bayesian estimation using the PPK model adding covariates of PPI and iron/magnesium oxide to the final model. As with the final model, best were the LSSs of C 6 in single sampling and of C 6 and C 8 in double sampling (Additional file 6 and Table 3). Taking into consideration the same results in Bayesian estimation and the fact that they are not always necessary for treating LN, we found it appropriate the decision to exclude concomitant PPI and iron /magnesium oxide from covariates.
The PPK model with 6 transit compartments in this study offered a good description of the absorption process. The model did not, however, improve with the addition of EHC, which causes the second peak in MPA concentration. In a previous study, EHC was reported to account for as much as 40% of MPA AUC 0-12 [4]. In this study, the second peak was seen in all but one patient, although peak values varied (Fig. 1A). Nevertheless, the EHC model did not significantly improve model fit. The failure of the EHC model may be because of the small number of patients, and also because of the lack of data on MPAG. MPA metabolizes into MPAG in the liver, and MPAG is thus involved in EHC [2]. In this study, however, MPAG concentrations were unknown and thus not taken into consideration. If the EHC model offered improved model fit, the prediction accuracy by MPA-BE would have been better.
As we mentioned earlier, AUC 0-12 of MPA is reported to be well associated with clinical outcomes [19,20] and is considered valuable for adjusting MMF dosage. As AUC 0-12 requires frequent blood sampling, possible surrogate markers such as MMF dose and trough concentration have been examined [31]. In this study, MMF dose did not correlate with MPA AUC 0-12 (r 2 = 0.14). Whereas MPA concentration at C 0 correlated with MPA AUC 0-12 (r 2 = 0.79) (Additional file 1), RMSPE was beyond the clinically acceptable range of 25% in evaluating the prediction accuracy of C 0 single blood sampling strategy by MLR. These results suggest that either MMF dose or C 0 MPA concentration may not work as surrogate markers for MPA AUC 0-12 . Further, no single sampling strategy by MLR showed RMSPE < 25%, suggesting that MPA concentration at any single time point is less likely usable as a surrogate marker. On the other hand, Bayesian estimation using the PPK model may allow a single blood sampling to predict AUC 0-12 with sufficient accuracy. In this study, the best 2 strategies of C 6 and C 8 by MAP-BE were good enough for clinical use (Table 3).
When double blood sampling is available, the prediction accuracy of LSS naturally increases. The best 4 strategies by both MAP-BE and MLR, except for double sampling at C 1 and C 3 , and double sampling at C 1 and C 4 by MAP-BE, were good enough for clinical use in this study. Although the best was the double sampling strategy of C 1 and C 4 by MLR, an LSS from MAP-BE may be preferable because of a few benefits. The MAP-BE approach might better estimate AUC 0-0.5 and C max , which this study demonstrated were associated with gastrointestinal involvement (Table 4), especially when only one or two blood samples are available. An LSS by MAP-BE would thus have a better potential of facilitating dose adjustment in consideration of both effectiveness and tolerability. In addition, estimation by MAP-BE can be flexible in terms of blood sampling timing, making LSSs by MAP-BE more practical from the clinical perspective [40]. The accuracy of MAP-BE estimation would increase if the study size is larger and the PPK model improves.
A plausible application of this study to clinical practice would be to determine the initial dosing using the PPK model and then adjust the dosage based on the Bayesian estimation using the measured values so that AUC 0-12 would be within the target range of 30-45 μg•h/mL. AUC 0-0.5 < 6.880 μg•h/mL and Cmax < 24.360 μg/mL would be desirable for minimizing gastrointestinal symptoms.
This study showed some limitations. First, this study did not fully examine all factors that may affect MPA pharmacokinetics. Prednisolone was administered in all patients and NSAIDs in only three patients: for these medications, comparisons were therefore impossible. Data on diet, which can reportedly affect MPA pharmacokinetics [47], were unavailable, as were data on pharmacogenetics. Secondly, this study did not take free MPA concentration into consideration, because the data were not available. In that free drug is therapeutically active, it would be desirable to evaluate free concentration in MPA. Measuring free drug concentration, though it may be still difficult in clinical practice, should be done in a future study. As another limitation, external validation of the PPK model was not conducted due to the small number of patients. Validation from a larger study is needed in the future. Especially, further examination is warranted in a study with a larger number of patients receiving PPI and iron/magnesium oxide. Finally, the study did not investigate the effectiveness of MMF, which should also be a subject for future study.

Conclusion
This study developed LSSs for Japanese LN patients receiving MMF so that the optimal MMF dosage could be individually determined based on estimated AUC 0-12 . When blood samples are taken twice, at 1 and 4 h after MMF administration, a strategy using the MLR approach would work best. When only a single blood sampling is allowed, however, a strategy using the MAP-BE approach with the PPK model developed in this study would be best with blood sampling at 6 h after administration. Although further studies are warranted, the LSS, whether single or double sampling, would facilitate safer, more effective, and individualized therapies for Japanese LN patients.